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1. Introduction 

Quantum key distribution (QKD) [H El [3] is the only quantum information protocol that 
found its way into practical applications, and is currently in a stage of early commercial 
development. There are two families of protocols that use fundamentally different 
resources. The original QKD protocol BB84 [4J and its variants transmit single photons 
(or approximations thereof), while the other family [5] perform measurements on pairs 
of entangled photons. A few years ago, entanglement-based QKD protocols were viewed 
as equivalent to BB84 [6] , and thus only of little interest for practical QKD due to their 
additional complexity. The new concept of device- independent QKD j7), and a returned 
awareness of classical side channels in prepare-and-send protocols revived interest in 
entanglement-based QKD schemes. Entangled photon pairs are efficiently prepared 
by spontaneous parametric down conversion (SPDC). Demonstrated for polarization- 
entangled pairs in 1995 [8], recent developments lead to the extremely bright sources 
available today (91110], so that entanglement-based QKD became a viable option. 

The first step in establishing a key in such a scheme is the assignment of 
photodetection events to entangled photon pairs. Due to their strong temporal 
correlation (down to a few 100 fs) in typical pair sources [H] , this assignment can be done 
via temporal coincidence identification. In typical laboratory experiments, as well as in 
early QKD implementations, a hardware channel was used to carry out this coincidence 
identification [T2l [T3j . Less hardware is required when coincidences are identified by 
comparing detection times given by good local clocks pH [15] or a central GPS time 
reference |16j . 

In this paper, we present an algorithm that relaxes the rather stringent reference 
clock quality requirements for such a coincidence identification so that conventional 
crystal oscillators can be used. In section [2j we outline the general problem and present 
a robust coincidence tracking scheme. Section [3] covers the algorithm to find an initial 
time offset as implemented in earlier experiments [TBI, ITT]. In sections HJ[6] we extend this 
scheme in the presence of a frequency difference between the clocks necessary to permit 
the use of clocks with lower accuracy. 

2. Photon pair identification with remote clocks 

The identification of pairs is straightforward in any context in which a hardware 
coincidence gate can be used; this is the case in laboratory-based experiments or field 
setups with a dedicated synchronization channel. 

The situation we address in this paper applies to cases where detection times of 
photons at the two distant locations [T5l [16j [18] are recorded, and coincidences are 
identified based on these time stamps (see figured]). This method requires stable and 
synchronous clocks used for the timestamping: A typical coincidence window r c is chosen 
to be slightly larger than the detector time jitter, which is on the order of Ins. The 
data acquisition for establishing a key out of measurements is supposed to run either 
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Figure 1. Setting of the problem. Detection times of photoevents from a correlated 
photon pair source and background are registered with respect to two local reference 
clocks at remote locations A, B. The true coincidences need then to be identified from 
the time sets {ti} and {t'j} on both sides. 



continuously, or at least for a few 100 seconds. To maintain two clocks synchronized 
within r c after a time of 100 s, a relative accuracy of 10 -11 is required, a specification 
that is met by commercial Rubidium clocks. For longer operation times, this still may 
be insufficient unless either a timing signal is transmitted on a separate channel, or the 
time reference is provided by a central source. 

Pair sources based on SPDC provide enough information in the streams of 
photodetection times {ti} and {t'j} that such accurate clocks should not be necessary. 
As long as the pair events are initially identified, the drift of the clocks can be tracked 
directly from the coincidence signal. For this to work reliably, the rate of pair events 
must be significantly larger than the one for accidental coincidences due to background 
photons in the same time window r c , which is also a necessary condition for a obtaining 
a secure key in QKD. . 

In its simplest form, a floating average of the time difference At = ti — t'j between 
true coincidence events can be used to track a drift of the reference time between the 
two sides. To illustrate this, and to evaluate the intrinsic clock stability necessary to 
follow the coincidence signature, we consider a realistic situation where the full width 
at half maximum of a coincidence time distribution due to detector jitter is = Ins. 
To estimate the center of this distribution with an uncertainty (one standard deviation) 
of 5t = 0.1 ns, we need to average time differences over about 

n=( -V + 1^19 (1) 

V2V2 In 25rJ 

coincidence events. Even for very low coincidence detection rates of 100 counts per 
second (cps), it takes less than 0.2 seconds to get a sufficient number of events. Over 
that period, the clock should not drift such that an event leaves the coincidence window, 
which translates into a relative frequency accuracy requirement of 10 -8 over 100 ms. 
More realistic coincidence detection rates of 1-10 keps require only a relative frequency 
accuracy of 10~ 7 to 10~ 6 over a period of 1 to 10 ms. Standard crystal oscillators easily 
exhibit a stability on that order, but may lack the accuracy. Thus, tracking the time 
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Figure 2. Effect of time offset and clock drift on photoevent sets. Trace (a) represents 
the event set {U} on side A, trace (b) an event set {t'j} on side B with a time offset 
AT, but the same reference clock frequency. Trace (c) illustrates a set {t'j} with an 
additional relative frequency difference Am between both reference clocks. 



difference in coincidences from a set of detection events permits to use these simpler 
reference oscillators during normal operation. 

Two problems are left for recovering the coincidences from time stamps derived 
with respect to two separate clocks: First, the detection instances at both sides will 
have an unknown time offset AT between them. This is mainly due to the absence 
of a common origin of time with a high enough resolution, and propagation over the 
physical distance between the two sides. As long as two reference clocks have the same 
frequency, AT can be found by looking at the cross correlation between the two timing 
signals. We will elaborate this in the next section. 

The second problem is related to the relative frequency difference between the two 
clocks due to a lack of accuracy. This is harder to solve, since the stream of time stamps 
{ti} and {t'j} on each side has no intrinsic time structure: Both signal and background 
events follow a Poisson distribution^]- 

The two problems of finding time- and frequency differences from coincidence signals 
in the presence of uncorrelated background events are illustrated in figure [2J Trace (a) 



| While it has been shown that the light emerging from SPDC processes should exhibit super-poissonian 
statistics [19], this is typically not observed in practical SPDC systems with photon counting detectors 
because it is washed out by multi-mode effects on a time scale way below the detector resolution. 
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shows a distribution of detection events {tj} on side A, trace (b) reflects the event stream 
on side B, assuming that there is only a time offset AT, but no frequency difference 
between the two reference clocks. Trace (c) shows an event stream in side B both under 
presence of a time offset and a frequency difference. For convenience, we describe the 
relative frequency difference by a quantity Au, such that the detection times t, t' on 
both sides due to identified photon pairs are connected via 

t' = (t + AT) • (f + Au). (2) 

We now estimate how accurately AT and Am need to be determined. In a practical 
QKD implementation, the two timestamping clocks are coarsely synchronized with 
conventional means (e.g. using an NTP protocol [20]), so it can be assumed that AT will 
not exceed a few 100 ms. A coincidence time window may be about 1 to 5 ns wide, fixing 
the uncertainty in AT to be small enough to start coincidence time tracking as sketched 
above. Thus, AT needs to be known with a precision of a few 10~ 9 , corresponding to an 
information of about 26 to 28 bit. For the tracking algorithm to take over, the relative 
frequency difference Au needs to be also known to an uncertainty of 10~ 8 to 10~ 6 . An 
upper bound for \Au\ can be chosen to match a typical accuracy of standard crystal 
oscillators (e.g. 10~ 4 ). Thus, Au of the two clocks needs to be found with a precision 
of 10~ 2 to 10 -4 , equivalent to an information of 7 to 14 bits. 



3. Finding the time offset 

We first explain the algorithm to find the time offset AT, assuming the two reference 
clocks run at the same frequency (Am = 0). Two streams of detection events {t{\ and 
{t'j} on both sides are translated into detection time functions 

a(t) = J2^(t-U), = (3) 

* 3 

The cross correlation between these two functions, 

c (r) = (a * 6)(r) := J a(t)b(t + r) dt , (4) 

has a peak at r = AT due to the correlated photodetection events on top of an 
unstructured but noisy base line from independent background detection events at both 
sides. The time offset AT is thus simply found by searching for the maximum in c(r) . In 
practice, c(r) is efficiently obtained from the timing sets via fast Fourier transformations 
(FFT) and their inverse, 

c{T) = r- 1 [F[a]-r[b]] , (5) 

with discrete arrays for a, b and c of length N (typically a power of 2). The high 
resolution necessary for AT (28 bits) renders a direct calculation impractical. It is 
possible, however, to obtain the coarse and fine part of AT separately with much smaller 
sample sizes. To illustrate how this works, we take the timing events {tj}, {t 1 -} captured 
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during an acquisition time T a , and map them onto the discrete arrays {a^}, with a 
time resolution 5t: 



y~ j^k,\jtil St) modNJ ; k = 0, . . . N — 1 , (6) 



and accordingly. This is an efficient process which requires visiting each entry 
t,i only once. The cross correlation array {c^} is obtained by the discrete version of 
equation ([5]), and its maximum located by a subsequent linear search in {c^}. If the cross 
correlation peak can be identified correctly, the result k max reflects AT up to a resolution 
5t, and modulo N5t. Thus, applying this method with two different resolutions 5t leads 
to a final AT with a resolution of 26 to 28 bit, while the individual FFTs are carried out 
at a moderate size of N = 2 19 or less. The complete code for this procedure is available 
as open source [2T] . 

It is beneficial to consider the influence of uncorrelated background events in this 
peak finding process. We assume a signal rate r s of true coincidences, and background 
rates r x and r 2 on both sides. The discrete arrays {a*,}, {bk} are built up from timestamps 
{ti}, {t'j} in a collection interval T a . The cross correlation peak will be made up 
by r s T a event pairs at the index k max , while the r^T^ background event pairs are 
homogeneously distributed over all N entries in {c^.} following a Poisson distribution. 
The peak can be identified with sufficient confidence if its statistical significance S, 
here defined as the ratio between the peak height above the base line and the standard 
deviation of the latter, 

8{k):= X^-= , (7) 
V (Cfe - c fc ) 2 

exceeds a certain numerical value. With the above rates, the peak value arising from 
signal pairs is 



xAwra/iv V r ^ 

If we approximate the fluctuations on the base line of {c^} by a Gaussian distribution, 
the probability e that a base line fluctuation gives rise to the largest value S max , and 
thus leads to a wrongly identified location of the cross correlation peak, is given by 

e = V (S max >S p )^^-(l- erf . (9) 

A numerical evaluation of this quantity (see table [T]) shows that for iV < 10 7 , S p > 6 
leaves less than 1% probability of misidentifying the peak. Since S max can be directly 
estimated out of {c^}, it forms a good basis to gauge the success of the peak finding 
procedure in practice. 

Care should be taken that events acquired over a time T a are uniformly distributed 
over the interval NSt in the binning procedure of equation (J6J). Specifically, T a /(N5t) 
should be an integer number. Otherwise, uncorrelated background events are subject 
to an effective envelope and do not lead to a flat base line in the cross correlation array, 
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Table 1. Connection between the probability e of wrong peak identification, bin 
number N and statistical significance S of a peak. 



e/N 


IO" 4 


IO" 5 


10~ 6 


lO" 7 


IO" 8 


IO" 9 


io- 10 


io- 11 


io- 12 


S 


3.72 


4.26 


4.75 


5.12 


5.61 


6.00 


6.36 


6.71 


7.03 




Figure 3. Cross correlation arrays {ck} of photo events acquired over T a w f.05s, 
normalized to a statistical significance S as defined in equation {7J with A*" = 2 19 . The 
kmax for two time resolutions St in (a) and (b) lead to a value AT = 53 599 f 60 ± 2 ns. 
All traces are sampled down by a factor 64. 

so determination of qT and subsequent peak finding becomes difficult. This problem can 
also be addressed by removing the lowest Fourier components in equation (jSJ) before the 
back transformation. 

From equation (jSJ) it can be seen that for given rates r 1; r 2 , and r s , the only 
way to increase the success probability is to increase the number of time bins N. 
For r\ = r 2 =100kcps, r s =lkcps, we need N « 360 000 < 2 19 to exceed S = 6. 
Furthermore, the frequency difference of the reference clocks needs to be very small. 
A time stretch T a Au between the two clocks over an acquisition time T a exceeding 
the targeted resolution St for AT reduces the statistical significance of the coincidence 
peak below a useful level. With parameters St=2ns and T a of a few seconds for the 
experiments carried out in [151 E] > reference clocks with a relative frequency difference 
Au < 10~ 9 were necessary, which were provided in the form of Rubidium oscillators. 

Figure [3] shows the result of typical correlation arrays {c^} (rescaled in terms of S) 
from an experiment with event rates r\ ~ 68 kcps, r 2 ~ 56 kcps and r s w 1280 cps. Here, 
N = 2 19 was chosen, and time resolutions St = 2 048 ns for the coarse, and St = 2 ns 
for the fine resolution. The peak exhibits S > 10 for both resolutions, and the resulting 
time offset is AT = 53 599 160 ± 2 ns. 
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4. Finding the time offset in the presence of a frequency difference 

The only reason to use reference clocks with a relative frequency accuracy better than 
10 -9 with the presented algorithms is to determine the initial time offset AT with a 
resolution on the order of 1 ns. Knowledge of the frequency difference to that accuracy 
and a reasonable stability is sufficient for tracking, so it is desirable to extract this 
information out of the timing events efficiently. 

Finding both the time and frequency difference from the recorded timing signals 
{ti}, {t'j} in the presence of uncorrelated background events is equivalent of identifying 
a line in a (ti,t' k ) plane of all possible pair events. This well-known pattern recognition 
problem is formally solved by the Hough transformation [22J , which maps the pair time 
distribution {(ti,t' k )} onto the parameter space {(AT, Au)}. As in the cross correlation 
method in the previous section, the pair (AT, Am) searched for is the peak coordinate 
in the parameter space. However, we did not find an equally efficient high-resolution 
solution as for the one-dimensional problem in section [3j 

A simpler method for determining Am is to estimate time offsets AT, AT 2 during 

relatively short acquisition intervals T a , shifted by a time T s > T a with the method 

described in the previous section. The change in time offsets between these probe 

intervals is connected with the relative frequency difference Am via 

ATi - AT 2 , . 

Au = — L- 1 . (10) 

J- s 

However, it is necessary to reliably obtain the time offsets on the two sampling intervals 
- which itself is only possible with clocks with a sufficiently small Am. We now evaluate 
under which conditions this cross correlation step will succeed in finding a time offset 
AT. 

For two clocks with Au = 0, the contribution of correlated events will all end up 
in a single time bin in the discrete correlation array {c^}. For Am ^ 0, the correlated 
events will spread out over roughly m = AuT a /St bin indices k. This reduces not only 
the statistical significance S for identifying a maximum, but also increases a timing 
uncertainty which in turn leads to an uncertainty in determining the frequency difference 
Am according to equation ffTUj) . 

In order to identify the correlation peak with sufficient confidence 1 — e according 
to equation ([9}, the statistical significance should exceed a threshold S t h ~ 6. For 
this, the timing resolution 5t may have to be increased, forcing the true coincidences 
in less bins k, up to to St = AuT a , or equivalently iV = 1/Aw. This, together with 
equation ([8]) for the statistical significance, leads to an expression for the maximally 
acceptable frequency difference 

Au max = q2 (11) 

for this strategy. In practice, the choice of a suitable size N for the correlation array 
{cj} has not to be done before the time-consuming step of the cross correlation in 
equation §5§. If with an initially chosen resolution St the peak is not found, the array 
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Figure 4. Correlation arrays for photoevents acquired with slightly detuned reference 
clocks. Traces (a) and (b) show arrays taken during acquisition time slots T a w 268 ms 
(Is for events at side B), separated by T s sa 1.074s and St — 2.048 /is - correlation 
peaks cannot be identified with sufficient significance. Traces (c) and (d) show the 
arrays after summing every 8 adjacent bins, revealing a moving correlation peak. All 
traces are sampled down. 



{cj} can be either re-partitioned in larger bins of width St', or equivalently exposed to a 
moving average procedure until a statistically significant correlation peak is identified. 
Once the time offset AT is known with an accuracy 5t', the relative frequency difference 
Au is known with an accuracy 

Si' T 

5u « V2%- = V2^r- • (12) 
T s NT S 

In the low signal limit, where N has to be chosen large enough to identify a peak at all, 
this results in a worst- case accuracy of Su = y/2T a /T s . 

We illustrate this method with experimental timing sets obtained from a non- 
optimized down conversion source with moderate background event rates (r\ ~ r 2 ~ 
77kcps, r s m 15kcps). Both timestamp units were referenced to crystal oscillators with 
a nominal frequency accuracy of 100 ppm. 

With two segments of T a = 2 28 ns ~ 268 ms and a separation of T s = 4T a m 1.074 s, 
convolution arrays {c^} were generated with a binning resolution of 5t = 2.048 fis. Since 
there were slowly varying changes in the background rates, the 20 lowest frequency 
entries (and their mirrors) were set to before back transformation to {c^}, resulting 
in a smooth base line. The results are shown in the top traces of figure HI without any 
significant correlation peaks. A subsequent re-binning with an effective width of St' = 
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16 384 ns reduced the noise level of the background sufficiently to allow the identification 
of the correlation peaks, resulting in time offsets of A7\ = AT = 374 505 ± 16 us, 
AT 2 = 374 292 ± 16 /is, and subsequently in Am = (1.98 ± 0.21) • 10~ 4 . 

5. Iterative procedure to decrease timing and frequency uncertainty 

The simple method for obtaining both AT and Au by analyzing correlations does 
typically not provide a sufficiently low uncertainty to start the tracking algorithm 
described in section [21 Therefore, additional steps are required. Knowledge of Au 
with uncertainty 5u < Au, the linear dependency of the timing uncertainty St from 
Au according to equation (Tl2l . and the time offset AT with some accuracy suggest an 
iterative method for this purpose: 

A set {t'j} is prepared from {t'j}, which is corrected with the initial values AT and 
Am (obtained as in section HI) via 

% = + AT) ■ (1 + Au) . (13) 

With this set and the original set {t{}, new values for Au and AT are obtained. The 
reduction in uncertainty 5u is given by the ratio T s /T a according to equation ffl2l . and 
is somewhat below an order of magnitude, or about 3 bit. This can be iterated, finally 
leading to values AT and Au with the targeted uncertainties (see appendix for the 
explicit algorithm). 

6. A faster algorithm for finding the fine time offset 

The iterative method in jointly finding AT and Am with sufficient accuracy converges 
only slowly because the time separation T s is typically not be very much larger than an 
acquisition time interval T a to keep the initial time for finding the coincidences low. 

Once initial values for Am and AT are found, an alternative algorithm can be 
used: We begin with event pairs sets {tj} and {t'j}, where the latter is corrected via 
equation ( TIBl similarly as before. If the timing resolution St for discretization is small 
enough (i.e., r li2 <5£ < 1), the arrays {a^} and {bj} are only sparsely populated. For 
t'j G [0,Xy with T b = 5t/5u, signal events lead to coincidences in time bins with the 
same time bin indices k = k! . The sparse population of the arrays {a^}, ensures 
then that the presence of a condition a^ = by=k = 1 is very likely due to a true 
coincidence. For those pair events, the instantaneous time difference At = ti — t'j due 
to an inaccurately known AT and Am can be determined with an accuracy limited only 
by the time resolution of the detection system (typically dominated by detector timing 
jitter). Analysis of instantaneous time differences At over the time interval T& finally 
reveals the parameters AT and Am with the intrinsic resolution of the system, after 
which the tracking algorithm can take over. 

A distribution of time differences generated with this method from the time stamps 
used in section H] is shown in figure [5]^a). With one more iteration of the correlation 
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Figure 5. (a) Time differences At between event pairs on both sides falling into the 
same time bin after pre-compensation with approximate AT, Am. A large fraction of 
the pairs appear on a line, with accidental coincidental pairs contributing to the noise 
of the figure. The differences fall in the range ±8t/2, and are known with a high 
precision, (b) Dropping adjacent pairs with excessive differences leaves a line which 
can be used to extract the final AT, Au. 

algorithm, values AT = 374 592.8 ±1.0 /xs and Au = 2.0113±0.0014- 10~ 4 were obtained 
to prepare the corrected set {t'A. The binning window for identifying coincidences was 
chosen as 5t = 1.024 /is. Out of the 2 19 bins for both sets {ti} and {t'j}, about 20 500 
were occupied with one entry, and 75 and 191 with two events each; the rest were empty, 
thus forming sparse arrays. 

One can visually identify a line structure, starting at about +0.25<5t, and increasing 
to about +0.55t towards the last of the 7839 coincidence candidates. Several of them 
are located away from this line, corresponding to bin pairs with accidental coincidences. 
The instantaneous time difference At of true coincidences increases only slowly with the 
binning index k, whereas the accidental coincidences can take arbitrary values. Thus, 
adjacent coincidence pairs with bin indices k < k! with a difference At^n— Ai^) in their 
instantaneous time difference exceeding a modulus of (tj(k')—ti(k))/N are likely to contain 
at least one accidental coincidence. In a cleaning step, such pairs of adjacent candidates 
are simply removed. This step left only a small number of 384 coincidence candidates 
in the list, apparently without any accidentals (see figure 15(b)). A linear fit with a 
model At = AT' + Au'ti with data from the remaining pairs returns offset correction 
parameters AT' = 260.5±0.05 ns and Am' = 4.270±0.002-10~ 7 . The error intervals from 
the fit appear overly optimistic, so we reconsider them, assuming a timing uncertainty 
of Ins in instantaneous measurements. We finally arrive at AT = 374 593 062 ± 1 ns and 
An = -200 789 ± 1.4- 10" 9 . 

These values are sufficient to start the tracking algorithm sketched in section [2j 
Figure [6] shows the evolution of the instantaneous time difference At, derived in a similar 
way as figure [5(b), but with the data corrected by the previously obtained constants 
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Figure 6. Time differences for coincidences after correction of the event times at side 
B with Au, AT. The variation can now be followed by a coincidence tracking scheme 
described in section [2] 

AT, Am for modeling the reference clock difference. One can recognize the slow drift of 
the reference oscillators, suggesting a stability around 10~ 8 on a time scale of a second. 

In conclusion, we presented an algorithm to remotely identify correlated photon 
pairs generated in a SPDC process from a stream of detection times without the need 
for a dedicated hardware channel, very stable and accurate pair of reference clocks or a 
central clock source, which may expose the classical infrastructure of a QKD system to 
a risk of compromise. This greatly reduces the technical complexity of entanglement- 
based QKD systems, making this part of an effort to simplify hardware by using intrinsic 
information in the photon pairs and bringing it closer to applications. 

This work is supported by the National Research Foundation & Ministry of 
Education, Singapore, and partly by a joint program of quantum information research 
between DSO and NUS. 

Appendix: Iterative algorithm for finding time- and frequency difference 

An explicit algorithm for obtaining time- and frequency differences AT an Au from 
time sets {tj}, {t'-} with high precision comprises the following steps: 

1. Choose the limits for the maximum expected frequency and time differences Au max 
and AT max ; the observed rates r 1; r 2 and an expected rate r s can be used to check 
with equation ffTTj) if this algorithm can be expected to be successful. 

2. Choose an acquisition time interval T a and a separation time interval T s , e.g. 
T s = 10T a . 
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3. Choose the smallest discretization time 5t that is compatible with a high chance 
of successfully identifying the correct peak in the cross correlation. Pair this 
with a suitable value of N for generating the arrays {aj} and {bk} according to 
equation (J6]). 

4. Generate the cross correlation array {c^} via FFT as in section [31 

5. Find the index k of the maximal value in {ck} and estimate its statistical significance 
S according to equation ([7]). 

6. If S is below a chosen significance limit S t h, half the size of the array {c^} by adding 
entries pairwise, and go back to step 5; this doubles the effective time resolution 
5t'. Otherwise, continue. 

7. Determine A7\ from the peak position k max and the effective time resolution St' of 
the last iteration of the previous step. 

8. With the time resolution 5t' in step 6, generate discrete arrays {a,-}, {bk} and {ck} 
for the second sampling interval, and determine from there AT 2 . 

9. Determine Am from A7\, AT 2 , and T s from equation ( flOl ). 

10. If 5t' in the last iteration is small enough to start the tracking as described in 
section El the algorithm is finished. 

11. Generate a modified set of event times {tj} according to t' = (t' + ATi) ■ (1 + Au) 

12. Choose the same N as in the last FFT, but reduce the time interval 5t by less than 
the expected gain in accuracy given by T s /T a /y/2; typically, this reduction factor 
would be 4 or 8 corresponding to 2 or 3 bits in accuracy gain. 

13. Generate new {a^}, {bj} and {ck} in the usual way from the original set {i,} and the 
modified set {t'j}; from the peak position in the new {c^}, determine the correction 
to ATi. Usually, this adds 2 or 3 bits in accuracy to ATi; proceed similarly for the 
correction to AT 2 . 

14. Continue with step 10. 
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